## Section 2 of the Online Appendix: including or not weekend days and holidays

library(dplyr)
library(ggplot2)
library(gridExtra)
library(grid)
library(reshape2)
library(grid)
library(rdrobust)
library(rdlocrand)
library(lubridate)


## Clear workspace

rm=list(ls())

##Opening the data:

load("ReplicationData.RData")


############
### 2010 ###
############

data2010$turnout <- as.numeric(data2010$VOTOU.1º.TURNO) - 1


#############################
## Discontinuity analysis ###
#############################

## 18-year olds

## Election Day: subsetting the 14-day window

data.18.rddED <- dplyr::filter(data2010, (dob >= "1992-09-27" & dob <= "1992-10-10") & turnout == 1)

DataRDDED18 <- data.18.rddED %>%
  filter(GRAU.INSTRUÇÃO != "Analfabeto") %>%
  group_by(dob) %>%
  summarise(voters = n(),
            turnout = unique(turnout))

DataRDDED18$daysToFrom <- seq(6, -7, by = -1)

##Using the rdrandinf package

rddED18 <- rdrandinf(DataRDDED18$voters, DataRDDED18$daysToFrom, wl = -7, wr = 6, seed = 50, ci = .05)

## % increase
rddED18$obs.stat/rddED18$sumstats[3,1]*100

## % increase: lower bound
rddED18$ci[1]/rddED18$sumstats[3,1]*100

## % increase: upper bound
rddED18$ci[2]/rddED18$sumstats[3,1]*100


##Identifying weekend days and holidays:

DataRDDED18$WeekendHolidays <- as.numeric(weekdays(DataRDDED18$dob) %in% c("Saturday", "Sunday"))

## Values for figure

meanTreat <- rddED18$sumstats[3,2]
meanControl <- rddED18$sumstats[3,1]

DataRDDED18$obliged <- as.numeric(DataRDDED18$daysToFrom >= 0)

lowerTreat <- meanTreat - qt(0.975,df=length(DataRDDED18$voters[DataRDDED18$obliged==1])-1)*rddED18$sumstats[4,2]/sqrt(length(DataRDDED18$voters[DataRDDED18$obliged==1]))

upperTreat <- meanTreat + qt(0.975,df=length(DataRDDED18$voters[DataRDDED18$obliged==1])-1)*rddED18$sumstats[4,2]/sqrt(length(DataRDDED18$voters[DataRDDED18$obliged==1]))

lowerControl <- meanControl - qt(0.975,df=length(DataRDDED18$voters[DataRDDED18$obliged==0])-1)*rddED18$sumstats[4,1]/sqrt(length(DataRDDED18$voters[DataRDDED18$obliged==0]))

upperControl <- meanControl + qt(0.975,df=length(DataRDDED18$voters[DataRDDED18$obliged==0])-1)*rddED18$sumstats[4,1]/sqrt(length(DataRDDED18$voters[DataRDDED18$obliged==0]))


##Graph:

graph <- ggplot(DataRDDED18,aes(x=dob,y=voters))+geom_point(aes(shape=factor(WeekendHolidays)))+geom_vline(xintercept=as.numeric(as.Date("1992/10/03")+0.5))+ylim(4200,6400)+ylab("Number of voters in 2010")+theme_bw()+theme(axis.title.x = element_text(vjust = 0),axis.title.y = element_text(vjust = 0.5))+scale_x_date(date_breaks = "2 day", date_labels = "%m-%d")+annotate('text',x=as.Date("1992/10/07"),y=5800,label="Not obliged to vote in 2010", hjust = 0.5)+annotate('text',x=as.Date("1992/09/30"),y=4600,label="Obliged to vote in 2010", hjust = 0.5)+theme(legend.position="none", axis.title.x=element_blank())+ggtitle("A. Including weekend days and holidays")

graph <- graph + geom_segment(aes(x=min(dob),xend = as.Date("1992/10/03"),y = meanTreat, yend=meanTreat)) + geom_segment(aes(x=max(dob),xend = as.Date("1992/10/04"), y = meanControl), yend=meanControl) + geom_segment(aes(x=min(dob),xend = as.Date("1992/10/03"),y = upperTreat, yend=upperTreat),color='grey')+ geom_segment(aes(x=min(dob),xend = as.Date("1992/10/03"),y = lowerTreat, yend=lowerTreat),color='grey')+ geom_segment(aes(x=max(dob),xend = as.Date("1992/10/04"),y = upperControl, yend=upperControl),color='grey') + geom_segment(aes(x=max(dob),xend = as.Date("1992/10/04"),y = lowerControl, yend=lowerControl),color='grey')


##Excluding weekend days:

DataRDDED18nowe <- DataRDDED18[ ! DataRDDED18$WeekendHolidays==1,]


##Using the rdrandinf package

rddED18nowe <- rdrandinf(DataRDDED18nowe$voters, DataRDDED18nowe$daysToFrom, wl = -7, wr = 6, seed = 50, ci = .05)

## % increase
rddED18nowe$obs.stat/rddED18nowe$sumstats[3,1]*100

## % increase: lower bound
rddED18nowe$ci[1]/rddED18nowe$sumstats[3,1]*100

## % increase: upper bound
rddED18nowe$ci[2]/rddED18nowe$sumstats[3,1]*100


## Values for figure

meanTreatnowe <- rddED18nowe$sumstats[3,2]
meanControlnowe <- rddED18nowe$sumstats[3,1]

lowerTreatnowe <- meanTreatnowe - qt(0.975,df=length(DataRDDED18nowe$voters[DataRDDED18nowe$obliged==1])-1)*rddED18nowe$sumstats[4,2]/sqrt(length(DataRDDED18nowe$voters[DataRDDED18nowe$obliged==1]))

upperTreatnowe <- meanTreatnowe + qt(0.975,df=length(DataRDDED18nowe$voters[DataRDDED18nowe$obliged==1])-1)*rddED18nowe$sumstats[4,2]/sqrt(length(DataRDDED18nowe$voters[DataRDDED18nowe$obliged==1]))

lowerControlnowe <- meanControlnowe - qt(0.975,df=length(DataRDDED18nowe$voters[DataRDDED18nowe$obliged==0])-1)*rddED18nowe$sumstats[4,1]/sqrt(length(DataRDDED18nowe$voters[DataRDDED18nowe$obliged==0]))

upperControlnowe <- meanControlnowe + qt(0.975,df=length(DataRDDED18nowe$voters[DataRDDED18nowe$obliged==0])-1)*rddED18nowe$sumstats[4,1]/sqrt(length(DataRDDED18nowe$voters[DataRDDED18nowe$obliged==0]))


##Graph:

graphnowe <- ggplot(DataRDDED18nowe,aes(x=dob,y=voters))+geom_point(aes(group=obliged))+geom_vline(xintercept=as.numeric(as.Date("1992/10/03")+0.5))+ylim(4200,6400)+xlab('Date of birth in 1992')+ylab("Number of voters in 2010")+theme_bw()+theme(axis.title.x = element_text(vjust = 0),axis.title.y = element_text(vjust = 0.5))+scale_x_date(date_breaks = "2 day", date_labels = "%m-%d")+annotate('text',x=as.Date("1992/10/07"),y=5800,label="Not obliged to vote in 2010", hjust = 0.5)+annotate('text',x=as.Date("1992/09/30"),y=4600,label="Obliged to vote in 2010", hjust = 0.5)+ggtitle("B. Excluding weekend days and holidays")

graphnowe <- graphnowe + geom_segment(aes(x=min(dob),xend = as.Date("1992/10/03"),y = meanTreatnowe, yend=meanTreatnowe)) + geom_segment(aes(x=as.Date("1992/10/10"),xend = as.Date("1992/10/04"), y = meanControlnowe), yend=meanControlnowe) + geom_segment(aes(x=min(dob),xend = as.Date("1992/10/03"),y = upperTreatnowe, yend=upperTreatnowe),color='grey')+ geom_segment(aes(x=min(dob),xend = as.Date("1992/10/03"),y = lowerTreatnowe, yend=lowerTreatnowe),color='grey')+ geom_segment(aes(x=as.Date("1992/10/10"),xend = as.Date("1992/10/04"),y = upperControlnowe, yend=upperControlnowe),color='grey') + geom_segment(aes(x=as.Date("1992/10/10"),xend = as.Date("1992/10/04"),y = lowerControlnowe, yend=lowerControlnowe),color='grey')

##Join both graphs:

graphcombined <- grid.arrange(graph, graphnowe, ncol=1)

ggsave(file="Online_Appendix_S2_plot_ED_2010.eps", plot=graphcombined)




## End-of-Year: subsetting the 14-day window

data.18.rddEoY <- dplyr::filter(data2010, (dob >= "1992-12-25" & dob <= "1993-01-07") & turnout == 1)

DataRDDEoY18 <- data.18.rddEoY %>%
  filter(GRAU.INSTRUÇÃO != "Analfabeto") %>%
  group_by(dob) %>%
  summarise(voters = n(),
            turnout = unique(turnout))

DataRDDEoY18$daysToFrom <- seq(6, -7, by = -1)

##Using the rdrandinf package

rddEoY18 <- rdrandinf(DataRDDEoY18$voters, DataRDDEoY18$daysToFrom, wl = -7, wr = 6, seed = 50, ci = .05)

## % increase
rddEoY18$obs.stat/rddEoY18$sumstats[3,1]*100

## % increase: lower bound
rddEoY18$ci[1]/rddEoY18$sumstats[3,1]*100

## % increase: upper bound
rddEoY18$ci[2]/rddEoY18$sumstats[3,1]*100


##Identifying weekend days and holidays:

DataRDDEoY18$WeekendHolidays <- as.numeric(weekdays(DataRDDEoY18$dob) %in% c("Saturday", "Sunday") | DataRDDEoY18$dob=="1992-12-24" | DataRDDEoY18$dob=="1992-12-25" | DataRDDEoY18$dob=="1992-12-31" | DataRDDEoY18$dob=="1993-01-01")

## Values for figure

meanTreat <- rddEoY18$sumstats[3,2]
meanControl <- rddEoY18$sumstats[3,1]

DataRDDEoY18$obliged <- as.numeric(DataRDDEoY18$daysToFrom >= 0)

lowerTreat <- meanTreat - qt(0.975,df=length(DataRDDEoY18$voters[DataRDDEoY18$obliged==1])-1)*rddEoY18$sumstats[4,2]/sqrt(length(DataRDDEoY18$voters[DataRDDEoY18$obliged==1]))

upperTreat <- meanTreat + qt(0.975,df=length(DataRDDEoY18$voters[DataRDDEoY18$obliged==1])-1)*rddEoY18$sumstats[4,2]/sqrt(length(DataRDDEoY18$voters[DataRDDEoY18$obliged==1]))

lowerControl <- meanControl - qt(0.975,df=length(DataRDDEoY18$voters[DataRDDEoY18$obliged==0])-1)*rddEoY18$sumstats[4,1]/sqrt(length(DataRDDEoY18$voters[DataRDDEoY18$obliged==0]))

upperControl <- meanControl + qt(0.975,df=length(DataRDDEoY18$voters[DataRDDEoY18$obliged==0])-1)*rddEoY18$sumstats[4,1]/sqrt(length(DataRDDEoY18$voters[DataRDDEoY18$obliged==0]))


##Graph:

graph <- ggplot(DataRDDEoY18,aes(x=dob,y=voters))+geom_point(aes(shape=factor(WeekendHolidays)))+geom_vline(xintercept=as.numeric(as.Date("1992/12/31")+0.5))+ylim(2700,4400)+ylab("Number of voters in 2010")+theme_bw()+theme(axis.title.x = element_text(vjust = 0),axis.title.y = element_text(vjust = 0.5))+scale_x_date(date_breaks = "2 day", date_labels = "%m-%d")+annotate('text',x=as.Date("1993/01/03"),y=2700,label="17-yo Cohort")+annotate('text',x=as.Date("1992/12/28"),y=2700,label="18-yo Cohort")+theme(legend.position="none", axis.title.x=element_blank())+ggtitle("A. Including weekend days and holidays")

graph <- graph + geom_segment(aes(x=min(dob),xend = as.Date("1992/12/31"),y = meanTreat, yend=meanTreat)) + geom_segment(aes(x=max(dob),xend = as.Date("1993/01/01"), y = meanControl), yend=meanControl) + geom_segment(aes(x=min(dob),xend = as.Date("1992/12/31"),y = upperTreat, yend=upperTreat),color='grey')+ geom_segment(aes(x=min(dob),xend = as.Date("1992/12/31"),y = lowerTreat, yend=lowerTreat),color='grey')+ geom_segment(aes(x=max(dob),xend = as.Date("1993/01/01"),y = upperControl, yend=upperControl),color='grey') + geom_segment(aes(x=max(dob),xend = as.Date("1993/01/01"),y = lowerControl, yend=lowerControl),color='grey')


##Excluding weekend days:

DataRDDEoY18nowe <- DataRDDEoY18[ ! DataRDDEoY18$WeekendHolidays==1,]


##Using the rdrandinf package

rddEoY18nowe <- rdrandinf(DataRDDEoY18nowe$voters, DataRDDEoY18nowe$daysToFrom, wl = -7, wr = 6, seed = 50, ci = .05)

## % increase
rddEoY18nowe$obs.stat/rddEoY18nowe$sumstats[3,1]*100

## % increase: lower bound
rddEoY18nowe$ci[1]/rddEoY18nowe$sumstats[3,1]*100

## % increase: upper bound
rddEoY18nowe$ci[2]/rddEoY18nowe$sumstats[3,1]*100

## Values for figure

meanTreatnowe <- rddEoY18nowe$sumstats[3,2]
meanControlnowe <- rddEoY18nowe$sumstats[3,1]

lowerTreatnowe <- meanTreatnowe - qt(0.975,df=length(DataRDDEoY18nowe$voters[DataRDDEoY18nowe$obliged==1])-1)*rddEoY18nowe$sumstats[4,2]/sqrt(length(DataRDDEoY18nowe$voters[DataRDDEoY18nowe$obliged==1]))

upperTreatnowe <- meanTreatnowe + qt(0.975,df=length(DataRDDEoY18nowe$voters[DataRDDEoY18nowe$obliged==1])-1)*rddEoY18nowe$sumstats[4,2]/sqrt(length(DataRDDEoY18nowe$voters[DataRDDEoY18nowe$obliged==1]))

lowerControlnowe <- meanControlnowe - qt(0.975,df=length(DataRDDEoY18nowe$voters[DataRDDEoY18nowe$obliged==0])-1)*rddEoY18nowe$sumstats[4,1]/sqrt(length(DataRDDEoY18nowe$voters[DataRDDEoY18nowe$obliged==0]))

upperControlnowe <- meanControlnowe + qt(0.975,df=length(DataRDDEoY18nowe$voters[DataRDDEoY18nowe$obliged==0])-1)*rddEoY18nowe$sumstats[4,1]/sqrt(length(DataRDDEoY18nowe$voters[DataRDDEoY18nowe$obliged==0]))


##Graph:

graphnowe <- ggplot(DataRDDEoY18nowe,aes(x=dob,y=voters))+geom_point(aes(group=obliged))+geom_vline(xintercept=as.numeric(as.Date("1992/12/31")+0.5))+ylim(2700,4400)+xlab('Date of birth in 1992/93')+ylab("Number of voters in 2010")+theme_bw()+theme(axis.title.x = element_text(vjust = 0),axis.title.y = element_text(vjust = 0.5))+scale_x_date(date_breaks = "2 day", date_labels = "%m-%d")+annotate('text',x=as.Date("1993/01/03"),y=2700,label="17-yo Cohort")+annotate('text',x=as.Date("1992/12/28"),y=2700,label="18-yo Cohort")+ggtitle("B. Excluding weekend days and holidays")

graphnowe <- graphnowe + geom_segment(aes(x=as.Date("1992/12/25"),xend = as.Date("1992/12/31"),y = meanTreatnowe, yend=meanTreatnowe)) + geom_segment(aes(x=as.Date("1993/01/07"),xend = as.Date("1993/01/01"), y = meanControlnowe), yend=meanControlnowe) + geom_segment(aes(x=as.Date("1992/12/25"),xend = as.Date("1992/12/31"),y = upperTreatnowe, yend=upperTreatnowe),color='grey')+ geom_segment(aes(x=as.Date("1992/12/25"),xend = as.Date("1992/12/31"),y = lowerTreatnowe, yend=lowerTreatnowe),color='grey')+ geom_segment(aes(x=as.Date("1993/01/07"),xend = as.Date("1993/01/01"),y = upperControlnowe, yend=upperControlnowe),color='grey') + geom_segment(aes(x=as.Date("1993/01/07"),xend = as.Date("1993/01/01"),y = lowerControlnowe, yend=lowerControlnowe),color='grey')


##Join both graphs:

graphcombined <- grid.arrange(graph, graphnowe, ncol=1)

ggsave(file="Online_Appendix_S2_plot_EoY_2010.eps", plot=graphcombined)

